home *** CD-ROM | disk | FTP | other *** search
/ SGI Hot Mix 17 / Hot Mix 17.iso / HM17_SGI / research / lib / primes.pro < prev    next >
Text File  |  1997-07-08  |  2KB  |  75 lines
  1. ;$Id: primes.pro,v 1.4 1997/01/15 03:11:50 ali Exp $
  2. ;
  3. ; Copyright (c) 1994-1997, Research Systems, Inc.  All rights reserved.
  4. ;       Unauthorized reproduction prohibited.
  5. ;+
  6. ; NAME:
  7. ;       PRIMES
  8. ;
  9. ; PURPOSE:
  10. ;       This function computes the first K prime numbers. The result is a 
  11. ;       K-element vector of type long integer.
  12. ;
  13. ; CATEGORY:
  14. ;       Statistics.
  15. ;
  16. ; CALLING SEQUENCE:
  17. ;       Result = Primes(K)
  18. ;
  19. ; INPUTS:
  20. ;       K:    A scalar of type integer or long integer that specifies the 
  21. ;             number of primes to be computed.
  22. ;
  23. ; EXAMPLE:
  24. ;       Compute the first 25 prime numbers.
  25. ;         result = primes(25)
  26. ;
  27. ;       The result should be:
  28. ;         [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31, 37, 41, 43, 47, $
  29. ;          53, 59, 61, 67, 71, 73, 79, 83, 89, 97]
  30. ;
  31. ; REFERENCE:
  32. ;       PROBABILITY and STATISTICS for ENGINEERS and SCIENTISTS (3rd edition)
  33. ;       Ronald E. Walpole & Raymond H. Myers
  34. ;       ISBN 0-02-424170-9
  35. ;
  36. ; MODIFICATION HISTORY:
  37. ;       Written by:  GGS, RSI, November 1994
  38. ;-
  39.  
  40. function primes, k
  41.  
  42.   ;Compute the first k prime numbers.
  43.  
  44.   prm = lonarr(k)
  45.  
  46.   prm[0] = 2L
  47.   n = 3L
  48.   count = 1L
  49.   prm[count] = 3L
  50.  
  51.   case2: count = count + 1L
  52.  
  53.   while(count lt k) do begin
  54.     case1:
  55.     n = n + 2L
  56.  
  57.     for ip = 1L, count do begin
  58.       q = n / prm[ip]
  59.       r = n mod prm[ip]
  60.       if r eq 0 then goto, case1  ;n is not prime.
  61.  
  62.       if q le prm[ip] then begin  ;n is prime.
  63.         prm[count] = n
  64.         goto, case2               ;compute next prime.
  65.       endif
  66.  
  67.     endfor 
  68.  
  69.   endwhile
  70.  
  71.   return, prm
  72.   prm = 0
  73.  
  74. end
  75.